This document contains code to fit the waggle dance model to each of our 20 different sites. The code here uses the wagglefit package, as well as some additional code stored in fit_data.R to help simplify things. Each site has its own section, with reused code. The final sections create map plots and summary statistic plots.
# set waggle dance duration in seconds as foraging distance for analysis
data <- read.csv("data/FullHBForagingData.csv")
alldata <- data %>%
select(date, site, duration.seconds) %>%
rename(foraging_distance = duration.seconds)
Provides a good fit on the data. All parameters look central and nicely covered.
target_site <- "BEL"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 1.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 0.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 8, -5.5, -1.8)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| BEL | collective | -180.8509 | 0.8917721 | 0.230719 | 5.150227 | 0.1709491 | 0.0348235 | 5 | 371.7018 | 0.0784314 | 0.904 |
| BEL | individual | -183.5011 | 1.0000000 | 4.495285 | NA | 0.0725148 | NA | 2 | 371.0021 | 0.0882353 | 0.812 |
all_sites[[target_site]]$fit
All parameters look central in the likelihood space and a nice fit is returned.
target_site <- "BFI"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 0.6)
bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.2, 3.5, -6.5, -2.)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| BFI | collective | -190.4988 | 0.04468 | 0.000001 | 2.76909 | 0.1643118 | 0.2488299 | 5 | 390.9975 | 0.0873362 | 0.329 |
| BFI | individual | -211.2967 | 1.00000 | 9.470271 | NA | 0.1000000 | NA | 2 | 426.5934 | 0.1834061 | 0.001 |
all_sites[[target_site]]$fit
The individual fit isnt great and the Bs paramater is increasing up to the boundry, indicating it can only really take on a straight line / exponential fit. The parameters for the collective model are fairly central in the likelihood space, the fit looks very good.
target_site <- "BLO"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-10, 5)
br_bnds <- c(1.0e-10, 5)
as_bnds <- c(1.0e-5, 5)
ar_bnds <- c(1.0e-10, 1.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 5505)
as_bnds <- c(1.0e-12, 0.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.2, 9.5, -8., -2.2)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| BLO | collective | -221.0260 | 0.2085372 | 1.785832 | 0.5709989 | 0.0677377 | 0.3591282 | 5 | 452.0519 | 0.0705882 | 0.775 |
| BLO | individual | -237.3888 | 1.0000000 | 5504.999895 | NA | 0.0001222 | NA | 2 | 478.7776 | 0.1117647 | 0.197 |
all_sites[[target_site]]$fit
Bs and Br go in oposite directions. E.g. Bs approaches 0 whilst Br approaches an every higher number.
Fit looks good.
target_site <- "BUR"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 500)
br_bnds <- c(1.0e-6, 500)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 1.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 0.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 4, -6.5, -2.2)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| BUR | collective | -118.3047 | 0.1121737 | 1e-06 | 500 | 0.1513434 | 0.0275123 | 5 | 246.6093 | 0.0670732 | 0.833 |
| BUR | individual | -143.7177 | 1.0000000 | 1e+01 | NA | 0.1000000 | NA | 2 | 291.4355 | 0.1768293 | 0.013 |
all_sites[[target_site]]$fit
The collective model roughly follows the individual model but is able to acheive a slightly improved fit to the shoulder, hence it provides a higher likelihood score. The AIC indicates this is overfitting, suggesting the individual model provides a more parsimonious explanation.
\(p\) does not approach 1 as one might expect, indicating there is overfitting and so a comparison with an individual model is required.
target_site <- "CAD"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 10)
ar_bnds <- c(1.0e-12, 10)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 2)
as_bnds <- c(1.0e-12, 0.8)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 1.8, -6, -1.2)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| CAD | collective | -48.27823 | 0.165521 | 1e-06 | 0.4246959 | 1.4564358 | 0.3861271 | 5 | 106.5565 | 0.0821918 | 0.943 |
| CAD | individual | -49.95312 | 1.000000 | 1e-06 | NA | 0.3854445 | NA | 2 | 103.9062 | 0.1095890 | 0.747 |
all_sites[[target_site]]$fit
The collective model provides the most parsimonious and best fit.
target_site <- "GIL"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 1.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 0.8)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 2.1, -8., -2)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| GIL | collective | -74.54094 | 0.26148 | 0.000001 | 1.8e-06 | 0.334516 | 0.7542132 | 5 | 159.0819 | 0.1302083 | 0.064 |
| GIL | individual | -102.30919 | 1.00000 | 2.036505 | NA | 0.345019 | NA | 2 | 208.6184 | 0.1979167 | 0.001 |
all_sites[[target_site]]$fit
The collective model provides the best fit to the data, but the proportion of scouts is high.
target_site <- "HER"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 100)
br_bnds <- c(1.0e-6, 100)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 1.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 2)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 4, -6.5, -1.5)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| HER | collective | -235.1801 | 0.6381663 | 0.0000010 | 1e-06 | 0.1777957 | 0.4291064 | 5 | 480.3603 | 0.1104651 | 0.230 |
| HER | individual | -242.3219 | 1.0000000 | 0.5484172 | NA | 0.1826291 | NA | 2 | 488.6438 | 0.1395349 | 0.072 |
all_sites[[target_site]]$fit
Collective model provides the best fit but falls under the tail and shoulder. The individual model strugles to find a good fit, probably due to the tail.
target_site <- "HHS"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 50)
br_bnds <- c(1.0e-6, 50)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 1.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 4, -6.5, -2.2)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| HHS | collective | -61.34797 | 0.0953653 | 0.000001 | 1e-06 | 0.1410585 | 0.5029664 | 5 | 132.6959 | 0.1333333 | 0.369 |
| HHS | individual | -83.74915 | 1.0000000 | 7.095173 | NA | 0.1251512 | NA | 2 | 171.4983 | 0.2444444 | 0.009 |
all_sites[[target_site]]$fit
Collective model provides a very good fit, whilst the individual model fails to find much traction. The proportion of scouts goes very low (~3%) suggesting the majority of the colony are following a small number of scouting individuals.
target_site <- "HOR"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 1.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 3, -9, -2.5)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| HOR | collective | -146.0403 | 0.0224734 | 0.000001 | 10 | 0.2126556 | 0.2025141 | 5 | 302.0806 | 0.0659898 | 0.335 |
| HOR | individual | -213.5839 | 1.0000000 | 4.342541 | NA | 0.2547302 | NA | 2 | 431.1677 | 0.1852792 | 0.000 |
all_sites[[target_site]]$fit
Again the collective model provides a good fit to the data, however, the individual model fits poorly, reducing to an exponential.
target_site <- "MAK"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 20)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# run individual model
individual_result <- fit_individual_model_to_data(data, bounds)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
# view individual model likelihood space to check bounds look ok
individual_result$llspace
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 4.2, -6.5, -2.5)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| MAK | collective | -84.51875 | 0.2093454 | 2.319366 | 0.0201563 | 0.1221729 | 0.4776196 | 5 | 179.0375 | 0.0909091 | 0.770 |
| MAK | individual | -98.02457 | 1.0000000 | 9.168369 | NA | 0.0961026 | NA | 2 | 200.0491 | 0.1919192 | 0.048 |
all_sites[[target_site]]$fit
Collective model provides the best fit, however it misses a large section of the shoulder for the tail. The \(bs\) parameter wants to go to zero, however when let go bellow 1e-6 the behaviour becomes very erratic and the fit deteriorates.
target_site <- "MEL"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 3., -7, -2.1)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| MEL | collective | -115.3164 | 0.0792579 | 1e-06 | 3e-06 | 0.2129400 | 0.3377206 | 5 | 240.6328 | 0.1626016 | 0.078 |
| MEL | individual | -135.7308 | 1.0000000 | 1e-06 | NA | 0.2520861 | NA | 2 | 275.4617 | 0.2439024 | 0.001 |
all_sites[[target_site]]$fit
Collective provides the best fit.
target_site <- "MPA"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 4., -7.2, -2.1)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| MPA | collective | -181.8729 | 0.566962 | 0.5783627 | 1e-06 | 0.1853641 | 0.419992 | 5 | 373.7458 | 0.0666667 | 0.888 |
| MPA | individual | -188.0378 | 1.000000 | 1.4738537 | NA | 0.1866656 | NA | 2 | 380.0756 | 0.1000000 | 0.429 |
all_sites[[target_site]]$fit
Colletive provides the best fit.
target_site <- "ROT"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 2., -5.5, -2.1)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| ROT | collective | -138.1294 | 0.3096621 | 1e-06 | 1e-06 | 0.3145301 | 0.4828992 | 5 | 286.2587 | 0.1494845 | 0.016 |
| ROT | individual | -153.3550 | 1.0000000 | 1e-06 | NA | 0.3537944 | NA | 2 | 310.7100 | 0.2010309 | 0.002 |
all_sites[[target_site]]$fit
Collective provides the best fit.
target_site <- "SAU"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 2.5, -6.5, -2.1)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| SAU | collective | -108.1848 | 0.3089541 | 0.000001 | 9.998126 | 0.2622360 | 0.1637328 | 5 | 226.3696 | 0.0606061 | 0.954 |
| SAU | individual | -115.5025 | 1.0000000 | 1.488528 | NA | 0.2710389 | NA | 2 | 235.0051 | 0.1287879 | 0.208 |
all_sites[[target_site]]$fit
Collective provides the best fit but misses a large section of the shoulder.
target_site <- "SOM"
# subset data for target site
data <- alldata %>%
filter(site == target_site) %>%
filter(foraging_distance < 6) # remove outlier distance
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 2., -5.5, -2.1)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| SOM | collective | -70.53597 | 0.0020792 | 1e+01 | 10 | 1.5000000 | 0.1486462 | 5 | 151.0719 | 0.1071429 | 0.528 |
| SOM | individual | -76.67133 | 1.0000000 | 1e-06 | NA | 0.3771972 | NA | 2 | 157.3427 | 0.1607143 | 0.080 |
all_sites[[target_site]]$fit
Collective provides the best fit.
target_site <- "SRA"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 3., -7.5, -2.1)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| SRA | collective | -123.3710 | 0.2468749 | 0.0960932 | 0.2081573 | 0.2126413 | 0.4919968 | 5 | 256.7419 | 0.0519481 | 0.979 |
| SRA | individual | -138.1142 | 1.0000000 | 3.0020146 | NA | 0.2133652 | NA | 2 | 280.2284 | 0.1428571 | 0.081 |
all_sites[[target_site]]$fit
Individual provides the best fit. Again, the proportion of scouts does not move towards 1. I.e. the collective model fails to reduce to the individual in the MLE fit.
target_site <- "STU"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 4, -6.5, -1.5)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| STU | collective | -155.6999 | 0.3518891 | 0.000001 | 1.028141 | 0.5214562 | 0.1557435 | 5 | 321.3997 | 0.0454545 | 1.000 |
| STU | individual | -156.5403 | 1.0000000 | 1.346564 | NA | 0.1612383 | NA | 2 | 317.0806 | 0.0545455 | 0.995 |
all_sites[[target_site]]$fit
Collective provides the best fit.
target_site <- "SWP"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 15)
br_bnds <- c(1.0e-6, 15)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 2., -6, -1.5)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| SWP | collective | -40.42382 | 0 | 0.1287557 | 1.274713 | 0.3653692 | 0.3841755 | 5 | 90.84765 | 0.1000000 | 0.722 |
| SWP | individual | -46.44077 | 1 | 0.0000010 | NA | 0.4236437 | NA | 2 | 96.88153 | 0.1444444 | 0.277 |
all_sites[[target_site]]$fit
Collective provides the best fit.
target_site <- "YAL"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.1)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 4., -8, -2.)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| YAL | collective | -256.2179 | 0.1526264 | 0.9940551 | 0.2527175 | 0.1890417 | 0.2865297 | 5 | 522.4357 | 0.0646552 | 0.715 |
| YAL | individual | -274.0674 | 1.0000000 | 0.5900943 | NA | 0.2131332 | NA | 2 | 552.1348 | 0.0948276 | 0.236 |
all_sites[[target_site]]$fit
Collective provides the best fit.
target_site <- "ZSL"
# subset data for target site
data <- alldata %>%
filter(site == target_site)
# set up bounds for the collective model
p_bnds <- c(0, 1.0)
bs_bnds <- c(1.0e-6, 10)
br_bnds <- c(1.0e-6, 10)
as_bnds <- c(1.0e-12, 1.5)
ar_bnds <- c(1.0e-12, 0.5)
collective_bounds <- rbind(
p_bnds, bs_bnds,
br_bnds, as_bnds,
ar_bnds
)
# set up bounds for the individual model
bs_bnds <- c(1.0e-6, 50)
as_bnds <- c(1.0e-12, 1.3)
individual_bounds <- rbind(
bs_bnds, as_bnds
)
# set coordinates for histogram subplot
subplot_coords <- c(0.1, 1.8, -8, -2)
all_sites[[target_site]] <- run_wagglefit_analysis(
target_site, data, collective_bounds, individual_bounds, subplot_coords
)
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
#> [1] "Itteration 1"
#> [1] "Itteration 2"
#> [1] "Itteration 3"
#> [1] "Itteration 4"
#> [1] "Itteration 5"
#> [1] "Itteration 6"
#> [1] "Itteration 7"
#> [1] "Itteration 8"
#> [1] "Itteration 9"
#> [1] "Itteration 10"
all_sites[[target_site]]$fit_result %>%
kbl() %>%
kable_classic(full_width = F)
| site | model | loglikelihood | p | bs | br | as | ar | k | AIC | ks_statistic | ks_pvalue |
|---|---|---|---|---|---|---|---|---|---|---|---|
| ZSL | collective | -150.2731 | 0.0522739 | 0.0000010 | 10 | 0.3736313 | 0.1997061 | 5 | 310.5463 | 0.0481013 | 0.710 |
| ZSL | individual | -181.7338 | 1.0000000 | 0.9650966 | NA | 0.4235067 | NA | 2 | 367.4675 | 0.1215190 | 0.003 |
all_sites[[target_site]]$fit
# save the analysis results for all sites
save(all_sites, file = "results/site_fit_results.Rdata")
# group all site results together
df <- map(all_sites, 1) %>%
bind_rows()
# save results
saveRDS(df, file = "results/site_fit_results.Rda")
# AIC plot
aic_plot <- df %>%
group_by(site) %>%
slice(which.min(AIC)) %>%
select(model) %>%
group_by(model) %>%
summarise(lowest_AIC = n()) %>%
ggplot(aes(x = model, y = lowest_AIC)) +
geom_bar(stat = "identity") +
labs(x = "Model", y = "Count") +
scale_y_continuous(breaks = seq(0, 20, by = 2)) +
theme(
text = element_text(size = 42)
)
ggsave(
plot = aic_plot,
filename = "results/figures/AIC_plot.png",
width = 90,
height = 110,
units = "mm",
dpi = 300
)
# ks plot
ks_plot_dist <- df %>%
ggplot(aes(x = ks_pvalue)) +
geom_histogram(bins = 10, binwidth = 0.1, col = "white") +
geom_vline(xintercept = 0.05, color = "red", linetype = "dashed") +
scale_y_continuous(breaks = seq(0, 12, by = 2)) +
labs(x = "KS P value") +
facet_wrap(~model, nrow = 3) +
theme(
text = element_text(size = 42),
strip.background = element_blank()
)
ggsave(
plot = ks_plot_dist,
filename = "results/figures/sites_ks.png",
width = 86,
height = 180,
units = "mm",
dpi = 300
)
ggsave(
plot = all_sites$STU$fit,
filename = "results/figures/STU.png",
width = 90,
height = 110,
units = "mm",
dpi = 300
)
ggsave(
plot = all_sites$ZSL$fit,
filename = "results/figures/ZSL.png",
width = 90,
height = 110,
units = "mm",
dpi = 300
)
Code to make the individual mal plots. These figures are created standalone but are combined in to facets manually in an image processor.
library(ggplot2)
library(ggrepel)
library(gridExtra)
library(ggsn)
library(sf)
library(rworldmap)
library(ggspatial)
library(rnaturalearth)
library(rnaturalearthdata)
full_data_path <- "data/FullHBForagingData.csv"
data_raw <- tibble(read.csv(full_data_path))
map_data <- data_raw %>%
select(
site, lat, lon
) %>%
distinct()
map_data <- df %>%
group_by(site) %>%
slice(which.min(AIC)) %>%
select(site, model) %>%
left_join(map_data, on = "site") %>%
mutate(col = ifelse(site %in% c("STU", "ZSL"), "1", "0"))
locations <- st_as_sf(
map_data,
coords = c("lon", "lat"), crs = 4326
)
# Extract selected sites for figure
selected_sites <- filter(map_data, site %in% c("STU", "ZSL")) %>%
mutate(
label = ifelse(site == "STU", "C (STU)", "D (ZSL)")
)
points_area <- st_bbox(locations)
worldmap <- ne_countries(scale = "large", returnclass = "sf")
# london area
london <- st_read(
"shapefiles/London_Ward.shp"
)
#> Reading layer `London_Ward' from data source
#> `/home/joe/Documents/PhD/wagglefit/analysis/shapefiles/London_Ward.shp'
#> using driver `ESRI Shapefile'
#> Simple feature collection with 649 features and 0 fields
#> Geometry type: POLYGON
#> Dimension: XY
#> Bounding box: xmin: 503568.2 ymin: 155850.8 xmax: 561957.5 ymax: 200933.9
#> Projected CRS: OSGB 1936 / British National Grid
london <- st_transform(
london,
CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0")
)
inset <- ggplot() +
geom_sf(
data = worldmap,
fill = "grey90",
color = "#4b4949d0"
) +
geom_sf(
data = london,
fill = "#ced0cffc",
lwd = 0
) +
coord_sf(
xlim = c(points_area[[1]] - 0.1, points_area[[3]] + 0.1),
ylim = c(points_area[[2]] - 0.1, points_area[[4]] + 0.1)
) +
geom_point(
data = map_data,
aes(x = lon, y = lat, shape = model, colour = model), size = 1.5
) +
geom_text(
data = selected_sites, aes(x = lon, y = lat, label = label),
nudge_x = c(0, -0.07), nudge_y = c(0.05, 0.05), size = 8
) +
annotation_north_arrow(
location = "tr", which_north = "true",
style = north_arrow_fancy_orienteering,
height = unit(15, "mm"),
width = unit(15, "mm"),
text_cex = 1.5
) +
annotation_scale(
location = "br",
text_cex = 1.5
) +
scale_shape_manual(values = c(1, 2)) +
scale_colour_manual(values = c("black", "red")) +
theme_nothing() +
scale_x_continuous(expand = c(0, 0)) +
scale_y_continuous(expand = c(0, 0)) +
labs(x = NULL, y = NULL) +
theme(
panel.border = element_rect(color = "black", fill = NA, size = .5)
)
inset
# make full plot
base <- ggplot(data = worldmap) +
geom_sf(
fill = "#c4cfc8",
color = "#4b4949d0",
lwd = 0.2
) +
coord_sf(
xlim = c(-11, 3),
ylim = c(49.5, 60)
) +
# geom_point(data = map_data, aes(x = lon, y = lat), size = 0.2) +
geom_rect(
aes(
xmin = points_area$xmin[[1]] - 0.1, xmax = points_area$xmax[[1]] + 0.1,
ymin = points_area$ymin[[1]] - 0.1, ymax = points_area$ymax[[1]] + 0.1
),
fill = NA,
colour = "black",
size = .02
) +
theme_nothing() +
scale_x_continuous(expand = c(0, 0)) +
scale_y_continuous(expand = c(0, 0)) +
labs(x = NULL, y = NULL) +
theme(
panel.border = element_rect(color = "black", fill = NA, size = .5)
)
base
merge_plot <- base +
annotation_custom(
ggplotGrob(inset),
xmin = 1,
xmax = 13,
ymin = 52.5,
ymax = 60
)
merge_plot
ggsave(
plot = merge_plot,
filename = "results/figures/site_map.png",
width = 90,
height = 110,
units = "mm",
dpi = 300
)
merge_plot_2 <- inset +
annotation_custom(
ggplotGrob(base),
xmin = -3.35,
xmax = 2,
ymin = 51.05,
ymax = 51.35
)
merge_plot_2
sites_model_plot <- plot_grid(
merge_plot_2, ks_plot_dist, all_sites$STU$fit, all_sites$ZSL$fit,
labels = c("A", "B", "C", "D"), label_size = 24
)
ggsave(
plot = sites_model_plot,
filename = "results/figures/sites_model_plot.png",
width = 183,
height = 190,
units = "mm",
dpi = 300
)
ggsave(
plot = sites_model_plot,
filename = "results/figures/sites_model_plot.svg",
width = 183,
height = 190,
units = "mm",
dpi = 300
)
ggsave(
plot = sites_model_plot,
filename = "results/figures/sites_model_plot.pdf",
width = 183,
height = 190,
units = "mm",
dpi = 300
)
model_sites_fits <- plot_grid(
merge_plot_2, ks_plot_dist,
labels = c("A", "B"), label_size = 24
)
ggsave(
plot = model_sites_fits,
filename = "results/figures/model_sites_fits.png",
width = 183,
height = 190,
units = "mm",
dpi = 300
)
stu_zsl_fit <- plot_grid(
all_sites$STU$fit, all_sites$ZSL$fit,
labels = c("A", "B"), label_size = 24
)
ggsave(
plot = stu_zsl_fit,
filename = "results/figures/stu_zsl_fit.svg",
width = 183,
height = 190,
units = "mm",
dpi = 100
)
# other map plots
inset <- ggplot() +
geom_sf(
data = worldmap,
fill = "grey90",
color = "#4b4949d0"
) +
geom_sf(
data = london,
fill = "#ced0cffc",
lwd = 0
) +
coord_sf(
xlim = c(points_area[[1]] - 0.1, points_area[[3]] + 0.1),
ylim = c(points_area[[2]] - 0.1, points_area[[4]] + 0.1)
) +
geom_point(
data = map_data,
aes(x = lon, y = lat, shape = model, colour = model), size = 1.5
) +
annotation_north_arrow(
location = "tr", which_north = "true",
style = north_arrow_fancy_orienteering,
height = unit(10, "mm"),
width = unit(10, "mm")
) +
annotation_scale(
location = "br",
text_cex = 3
) +
scale_shape_manual(values = c(1, 2)) +
scale_colour_manual(values = c("black", "red")) +
theme_nothing() +
scale_x_continuous(expand = c(0, 0)) +
scale_y_continuous(expand = c(0, 0)) +
labs(x = NULL, y = NULL) +
theme(
panel.border = element_rect(color = "black", fill = NA, size = .5),
text = element_text(size = 42)
)
inset
# make full plot
base <- ggplot(data = worldmap) +
geom_sf(
fill = "#c4cfc8",
color = "#4b4949d0",
lwd = 0.2
) +
coord_sf(
xlim = c(-11, 15),
ylim = c(49.5, 60)
) +
geom_point(data = map_data, aes(x = lon, y = lat), size = 0.2) +
geom_rect(
aes(
xmin = points_area$xmin[[1]] - 0.1, xmax = points_area$xmax[[1]] + 0.1,
ymin = points_area$ymin[[1]] - 0.1, ymax = points_area$ymax[[1]] + 0.1
),
fill = NA,
colour = "black",
size = .02
) +
theme_nothing() +
scale_x_continuous(expand = c(0, 0)) +
scale_y_continuous(expand = c(0, 0)) +
labs(x = NULL, y = NULL) +
theme(
panel.border = element_rect(color = "black", fill = NA, size = .5)
)
merge_plot <- base +
annotation_custom(
ggplotGrob(inset),
xmin = 1,
xmax = 13,
ymin = 52.5,
ymax = 60
)
merge_plot
ggsave(
plot = merge_plot,
filename = "results/figures/site_map.png",
width = 183,
height = 190,
units = "mm",
dpi = 300
)
sites_model_plot <- plot_grid(
all_sites$STU$fit, all_sites$ZSL$fit, aic_plot, ks_plot_dist,
labels = c("A", "B", "C", "D"), label_size = 42
)
sites_model_plot
ggsave(
plot = sites_model_plot,
filename = "results/figures/results_model_plot.png",
width = 183,
height = 190,
units = "mm",
dpi = 300
)